import numpy as np
import random
from scipy import sparse
import matplotlib.pyplot as plt
from ctypes import *
import numpy.ctypeslib as npct
import ctypes
from collections import Counter
import networkx as ne
import seaborn as sns
class society_network:
    def __init__(self,a,b,N,m,p,Np,Nf):
        self.touch_d = []
        self.N = N  # 模拟总人数
        self.a = a  # 最小家庭人口数
        self.b = b  # 最大家庭人口数
        self.touch_sum = np.zeros((N, N))
        self.touch_record = np.zeros((N, N)) - 1  # touch_record现在用来存储最近x天内有过接触行为,且初始值置为-1
        self.family_num = []  # 每一类家庭人口的数目
        self.family_mem = []  # 家庭信息列表,记录了每个家庭的首位个体序号和其家庭人数
        self.public_num = 100  # 公共网络抽样数量
        self.random_list = random.sample(range(self.N), self.public_num)
        self.empty_list = []  # 当前空位序号
        self.sample_fam_p = 0.001  # 初始化空位所占比例
        self.isolated_people = []  # 当前被隔离者
        # 小世界网络参数
        self.m = m
        self.p = p
        # 网络活跃时长
        self.Np = Np
        self.nf = Nf
        # 网络平均度
        self.mean_quantities = []

    def family_pro(self):
        '''
        划分家庭信息，并存入family_mem,family_num
        并对非整除情况，对其余数进行处理
        '''
        a, b, N = self.a, self.b, self.N
        sum = (a + b) * (b - a + 1) / 2
        org_num = int(N / sum)
        last_num = self.N - (org_num * sum)
        for i in range(a, b + 1):
            self.family_num.append(org_num)
            if last_num and last_num >= i:
                last_num = last_num - i
                self.family_num[i - a] = self.family_num[i - a] + 1
            # 判断此时余数是否大于当前i，若小于则此时的余数为一个家庭
            elif last_num:
                self.family_num[int(last_num - a)] = self.family_num[int(last_num - a)] + 1
                last_num = 0
        # 将家庭信息存储，存储每个家庭的头号和人口数
        i = 0
        for s in range(self.a, self.b + 1):
            for n in range(self.family_num[s - self.a]):
                self.family_mem.append([i, s])
                i = i + s

    def small_world_NW(self, N, m, p):
        '''
        使用ctypes调用小世界网络c程序
        :param N: 小世界网络个体数
        :param m: m
        :param p: 重连概率
        :return: 返回其接触矩阵
        # '''
        # FX = ctypes.cdll.LoadLibrary
        # lib_ctype = FX("./ketizu.dll")
        # reg = np.zeros((N, N), dtype=np.int32)
        # lib_ctype.small_world.argtypes = [npct.ndpointer(dtype=np.int, shape=(N, N)), c_int, c_int, c_double]
        # lib_ctype.small_world(reg, c_int(N), c_int(m), c_double(p))
        reg = ne.watts_strogatz_graph(self.N, self.m, self.p)
        reg = np.array(ne.adjacency_matrix(reg).todense())
        return reg

    def public_net(self):
        '''
        公共网络接触模拟：
        将个体分别投入到公共网络中，公共网络包含四个子网络，均用小世界网络算法模拟接触，
        将四个子网络的接触行为汇总到矩阵touch_sum
        '''
        if self.public_num > 0:
            # 公共网络划分
            size = self.public_num // 4
            list_a = self.random_list[:size]
            list_b = self.random_list[size:2 * size]
            list_c = self.random_list[2 * size:3 * size]
            list_d = self.random_list[3 * size:self.public_num]
            # sample个体投入公共网络
            TA = self.small_world_NW(size, self.m, self.p)
            TB = self.small_world_NW(size, self.m, self.p)
            TC = self.small_world_NW(size, self.m, self.p)
            TD = self.small_world_NW(size, self.m, self.p)
            for i in range(size):
                for j in range(size):
                    if TA[i, j]:
                        p = list_a[i]
                        q = list_a[j]
                        self.touch_sum[p, q] = 1
                        self.touch_sum[q, p] = 1
                    if TB[i, j]:
                        p = list_b[i]
                        q = list_b[j]
                        self.touch_sum[p, q] = 1
                        self.touch_sum[q, p] = 1
                    if TC[i, j]:
                        p = list_c[i]
                        q = list_c[j]
                        self.touch_sum[p, q] = 1
                        self.touch_sum[q, p] = 1
                    if TD[i, j]:
                        p = list_d[i]
                        q = list_d[j]
                        self.touch_sum[p, q] = 1
                        self.touch_sum[q, p] = 1

    def family_net(self):
        '''
        家庭网络的接触模拟：
        按照family_mem,family_num中存储的家庭信息，将家人之间设置为完全接触
        '''
        for item in self.family_mem:
            for q in range(item[0], item[0] + item[1] - 1):
                for m in range(q + 1, item[0] + item[1]):
                    self.touch_sum[m][q] = 1
                    self.touch_sum[q][m] = 1

    def people_in(self, i, s):
        '''
        迁入者模拟：
        在当前空位中抽取一定序号作为迁入者，并从中删除这些个体序号，将这些序号加入到
        enter_list中，并对这些序号进行车厢模拟接触。
        :param i: 迁入人数
        :param s: 车厢接触半径
        :return: 本次的迁入者序号，由主程序并赋予其状态属性
        '''
        enter_list = []
        if i > len(self.empty_list):
            print("当前迁入人数超出空位数！")
        else:
            # 从empty——list中抽取迁入者序号，并删除
            r_list = random.sample(range(len(self.empty_list)), i)
            for r in r_list:
                enter_list.append(self.empty_list[r])
                # 将迁入个体初始置为0
                self.touch_sum[r][:] = 0
                self.touch_sum[:][r] = 0
            # 删除元素
            self.empty_list = [self.empty_list[i] for i in range(0, len(self.empty_list), 1) if i not in r_list]
            # 迁入者的车厢接触模拟
            for a in range(0, len(enter_list), s):
                for i in range(s):
                    if a + i < len(enter_list):
                        self.touch_sum[enter_list[a]][enter_list[a + i]] = 1
            return enter_list

    def people_out(self, o):
        '''
        抽取一定数量的迁出者序号并加入self.empty_list(从当前非空位序号中抽取），将空位置为-1；
        :param o: 迁出人数
        '''
        remain_list = list(set(range(self.N)).difference(set(self.empty_list)))
        o_list = random.sample(range(len(remain_list)), o)
        for o in o_list:
            self.empty_list.append(remain_list[o])
        # 迁出人员置为-1
        for i in self.empty_list:
            self.touch_sum[i, :] = -1
            self.touch_sum[:, i] = -1

    def touch_record_update(self):
        '''
        用touch_record矩阵现在用来存储最近一次接触行为在x天前
        对该矩阵处理若本次无接触则该矩阵对应值+1，若有接触行为则置为1；
        '''
        for i in range(self.N):
            for j in range(self.N):
                if self.touch_sum[i][j] == 0 and i != j:
                    if self.touch_record[i][j] > 0:
                        self.touch_record[i][j] += 1
                    else:
                        pass
                if self.touch_sum[i][j] == 1:
                    # 如果当前个体是被隔离者，则当作本次无接触处理，即忽略接触行为
                    if (i in self.isolated_people) or (j in self.isolated_people):
                        self.touch_record[i][j] += 1

                    else:
                        self.touch_record[i][j] = 1
                if self.touch_record[i][j] == 15:
                    self.touch_record[i][j] = 0
    def Network_initialization(self):
        '''进行家庭和空位初始化'''
        self.family_pro()
        #空位初始化
        # for item in self.family_mem:
        #     p = np.random.uniform(0, 1, 1)
        #     if p < self.sample_fam_p:
        #         for i in range(item[1]):
        #             self.empty_list.append(item[0] + i)
        #     else:
        #         pass
        # for i in self.empty_list:
        #     self.touch_sum[i, :] = -1
        #     self.touch_sum[:, i] = -1
        return self.touch_sum

    def Network_activity(self):
        '''
        按一定的网络活跃时间分别调用两种接触网络运行
        :return:
        '''
        for i in range(self.Np):
            self.public_net()
        for f in range(self.nf):
            self.family_net()
        return self.touch_sum

    def step(self, o, i, s):
        '''
        迭代函数：
        :param o: 本次迁出人数
        :param i: 本次迁入人数
        :param s: 车厢感染半径
        :return:
        '''
        self.touch_sum = np.zeros((self.N, self.N))
        self.Network_activity()
        # self.people_out(o)
        # self.people_in(i, s)
        self.touch_record_update()

    def update_everyday(self, isolation_people):
        self.isolated_people += isolation_people
        self.touch_sum = np.zeros((self.N, self.N))
        self.Network_activity()
        self.touch_record_update()
        return self.touch_sum

    def avg_quantities(self):
        sum_p = np.zeros((self.N,))
        matrix = self.touch_record
        for i in range(self.N):
            for j in range(self.N):
                if matrix[i, j] > 0:
                    sum_p[i] += 1
                else:
                    sum_p[i] += 0

        avg_t = np.sum(sum_p) / len(sum_p)
        self.touch_d.append(avg_t)

    def dufenbu(self, matrix):  # 测试当前网络度分布
        '''
        测量当前矩阵度分布
        :param matrix: 输入矩阵
        :return: 分别返回度的概率，度的值和度的频次列表
        '''
        quantities = np.sum(matrix, axis=1)
        pindu = Counter(quantities)
        pindu = dict(sorted(pindu.items()))
        p_num = list(pindu.values())
        p_k = [x / sum(p_num) for x in p_num]
        x_label = list(pindu.keys())
        return p_k,x_label,quantities
    def juleixi(self,G, i):
        '''

        :param i: 节点序号
        :return: 该节点的聚类系数
        '''
        (m, n) = np.shape(G)
        neighbor_list = []  # i的邻居序号
        c = 0  # i的邻居对为直连的个数
        for j in range(m):
            if G[i, j] == 1:
                neighbor_list.append(j)
        for neighbor1 in neighbor_list:
            for neighbor2 in neighbor_list:
                if G[neighbor1, neighbor2] == 1:
                    c = c + 1
        node_sum = len(neighbor_list) * (len(neighbor_list) - 1)
        xishu = c / node_sum  # 聚类系数
        return xishu
if __name__ == "__main__":
    def tiaocan(args,args_list):
        N=1000
        p=0.6
        m=10
        np=1
        nf=1
        if args=='m':
            society1 = society_network(1,5,N,args_list[0],p,np,nf)#a,b,N,m,p,np,nf
            society3 = society_network(1,5,N,args_list[1],p,np,nf)
            society5  = society_network(1,5,N,args_list[2],p,np,nf)
            society7 = society_network(1,5,N,args_list[3],p,np,nf)
        elif args=='p':
            society1 = society_network(1, 5, N, m, args_list[0], np, nf)  # a,b,N,m,p,np,nf
            society3 = society_network(1, 5, N, m, args_list[1],np, nf)
            society5 = society_network(1, 5, N, m, args_list[2],np, nf)
            society7 = society_network(1, 5, N, m, args_list[3],np, nf)
        elif args=='np':
            society1 = society_network(1, 5, N, m, p, args_list[0], nf)  # a,b,N,m,p,np,nf
            society3 = society_network(1, 5, N, m, p, args_list[1], nf)
            society5 = society_network(1, 5, N, m, p, args_list[2], nf)
            society7 = society_network(1, 5, N, m, p, args_list[3], nf)
        society1.Network_initialization()
        society3.Network_initialization()
        society5.Network_initialization()
        society7.Network_initialization()
        for day in range(15):
            society1.step(1, 1, 10)
            society1.avg_quantities()
            society3.step(1, 1, 10)
            society3.avg_quantities()
            society5.step(1, 1, 10)
            society5.avg_quantities()
            society7.step(1,1,10)
            society7.avg_quantities()
            print(',',end='')
        paint_avgquantities(args,args_list,society1,society3,society5,society7)
    #画图1
    def paint_avgquantities(args,args_list,s1,s3,s5,s7):
        plt.plot(range(15),s1.touch_d,color='red',label='{}={}'.format(args,args_list[0]))
        plt.plot(range(15),s3.touch_d,color='green',label='{}={}'.format(args,args_list[1]))
        plt.plot(range(15), s5.touch_d, color='blue',label='{}={}'.format(args,args_list[2]))
        plt.plot(range(15), s7.touch_d, color='pink', label='{}={}'.format(args,args_list[3]))
        plt.rcParams['font.sans-serif'] = ['SimHei']
        plt.rcParams['axes.unicode_minus'] = False
        plt.title("np={},p={},nf={}".format(society1.np,society1.p,society1.nf))
        plt.xlabel("天数")
        plt.ylabel("平均每个人的有效接触人数")
        plt.legend()
        plt.show()
    tiaocan('m',[6,8,10,12])

    # pk_1,x1,q1=society1.dufenbu(society1.touch_sum)
    # pk_2,x2,q2=society3.dufenbu(society3.touch_sum)
    # pk_3,x3,q3=society5.dufenbu(society5.touch_sum)
    # pk_4,x4,q4=society7.dufenbu(society7.touch_sum)
    # sns.kdeplot(q1, shade=True, color='r',label='p=0.01')
    # sns.kdeplot(q2, shade=True, color='blue',label='p=0.1')
    # sns.kdeplot(q3, shade=True, color='green',label='p=0.4')
    # sns.kdeplot(q4, shade=True, color='pink',label='p=0.8')
    # plt.xticks(range(25))
    # plt.yticks([z/10 for z in range(5)])
    # plt.legend()
    # plt.show()
